Click here to view/hide the code
## data viz for humidity yerawise
pacman::p_load("ggcorrplot","directlabels","pglm", "scales")
suppressPackageStartupMessages(library(tidyverse))Humidity-CVD
Data loading
Packages
Correlation Analysis
Click here to view/hide the code
# Correlation analysis
combined_data <- combined_data |>
rename(maximum_25=high_25,
maximum_975=high_975,
minimum_25=low_25,
minimum_975=low_975,
average_25=avg_25,
average_975=avg_975,
median_25=med_25,
median_975=med_975,)
correlation_matrix <- cor(combined_data |> select(average_25, average_975, median_25,median_975,mode_25,mode_975,
maximum_25,maximum_975,
minimum_25,minimum_975,variation_25,variation_975), method = "pearson")
#correlation_matrix
# visualize as heatmap
# Compute a matrix of correlation p-values
p.mat <- cor_pmat(correlation_matrix)Click here to view/hide the code
ggcorrplot(correlation_matrix,
outline.color = "white",
type = "lower",
ggtheme = ggplot2::theme_gray,
colors = c("#6D9EC1", "white", "#E46726"),
lab = TRUE,
p.mat = p.mat,
insig = "blank")
This matrix shows pairwise correlations among various humidity metrics (average, median, mode, minimum, maximum, variation) at the 2.5th and 97.5th percentiles. Only statistically significant correlations are displayed, indicating strong associations among several metrics. Based on this analysis, both 2.5 and 97.5 percentiles of the average humidity, and 2.5th percentile of the variation in humidity were selected for further analysis. This is to minimize multicollinearity and retain representative variability in the downstream panel regression.
Click here to view/hide the code
# averages
combined_data |>
ggplot(aes(x = factor(Year), group=1)) +
geom_line(aes(y=average_25),linetype="dashed", color="midnightblue")+
geom_line(aes(y=average_975),linetype="dashed", color="midnightblue")+
geom_ribbon(aes(ymin=average_25, ymax=average_975),fill="midnightblue", alpha=0.4)+
facet_wrap(~State)+
theme_minimal(base_size = 20)+
labs(y="humidity (%)", x="Year",
title = "Average humidity",
subtitle = "2.5th and 97.5th percentiles")+
theme(axis.text.x = element_text(angle = 90, hjust = 1, size = 12))
# variations
combined_data |>
ggplot(aes(x = factor(Year), group=1)) +
geom_line(aes(y=variation_25),linetype="dashed", color="seagreen")+
geom_line(aes(y=variation_975),linetype="dashed", color="seagreen")+
geom_ribbon(aes(ymin=variation_25, ymax=variation_975),fill="seagreen", alpha=0.4)+
#geom_ribbon(aes(ymin=0, ymax=variation_25),fill="seagreen", alpha=0.4)+
facet_wrap(~State)+
theme_minimal(base_size = 20)+
labs(y="humidity (%)", x="Year",
title = "Humidity Variations",
subtitle = "2.5th and 97.5th percentiles")+
theme(axis.text.x = element_text(angle = 90, hjust = 1, size = 12))
Figure 2 (a) illustrates the trends in the 2.5th and 97.5th percentiles of daily humidity for each year across Indian states from 2011 to 2019. The 2.5th percentile of annual daily humidity remains consistently below 25% across all states. The lowest value was observed in Gujarat in 2011 (13.8%), while the highest value at this percentile was recorded in Tamil Nadu in 2019 (22.1%). This indicates a limited range of variation at the lower end of the humidity spectrum across states. The 97.5th percentile of annual daily humidity showed a wider range of variability across states, with the lowest value recorded in Odisha in 2011 (37.3%) and the highest value in Uttarakhand in 2019 (87.2%). While most states consistently recorded values below 75%, certain states crossed this threshold in specific years. For instance, Uttarakhand, Sikkim, Rajasthan, Punjab, Haryana, Himachal Pradesh, Delhi, Bihar, and Arunachal Pradesh recorded 97.5th percentile values above 75% in 2019. Similarly, northeastern states like Nagaland, Meghalaya, and Manipur crossed this threshold in 2013. These observations indicate episodic peaks in maximum humidity levels in certain states, particularly during the later years of the study period.
Figure 2 (b) illustrates the trends in the 2.5th and 97.5th percentiles of annual variations in daily humidity across Indian states from 2011 to 2019. Variations at the 2.5th percentile ranged from 4.6% in Delhi in 2016 to 10% in Andhra Pradesh in 2016, indicating limited variability at the lower bound across states and years. Most states exhibited stable patterns at the 2.5th percentile, with values generally below 10%. In contrast, the 97.5th percentile of annual humidity variations displayed more pronounced inter-state and temporal differences, ranging from 18.6% in Rajasthan in 2017 to 49.7% in Himachal Pradesh in 2019.
Click here to view/hide the code
model_CVD_Death_avg25 <- pglm(CVD_Death ~ average_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Death_avg25)
confint(model_CVD_Death_avg25)[, 1]
confint(model_CVD_Death_avg25)[, 2]
irr_avg25_death <- exp(coef(model_CVD_Death_avg25))
irr_avg25_death
irr_lci_avg25_death <- exp(confint(model_CVD_Death_avg25)[, 1])
irr_uci_avg25_death <- exp(confint(model_CVD_Death_avg25)[, 2])
irr_lci_avg25_death
irr_uci_avg25_death
model_CVD_Death_avg975 <- pglm(CVD_Death ~ average_975,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Death_avg975)
confint(model_CVD_Death_avg975)[, 1]
confint(model_CVD_Death_avg975)[, 2]
irr_avg975_death <- exp(coef(model_CVD_Death_avg975))
irr_avg975_death
irr_lci_avg975_death <- exp(confint(model_CVD_Death_avg975)[, 1])
irr_uci_avg975_death <- exp(confint(model_CVD_Death_avg975)[, 2])
irr_lci_avg975_death
irr_uci_avg975_death
model_CVD_Death_var25 <- pglm(CVD_Death ~ variation_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Death_var25)
irr_var25_death <- exp(coef(model_CVD_Death_var25))
irr_var25_death
irr_lci_var25_death <- exp(confint(model_CVD_Death_var25)[, 1])
irr_uci_var25_death <- exp(confint(model_CVD_Death_var25)[, 2])
irr_lci_var25_death
irr_uci_var25_death
### CVD_Incidence
model_CVD_Incidence_avg25 <- pglm(CVD_Incidence ~ average_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Incidence_avg25)
irr_avg25_incidence <- exp(coef(model_CVD_Incidence_avg25))
irr_avg25_incidence
irr_lci_avg25_incidence <- exp(confint(model_CVD_Incidence_avg25)[, 1])
irr_uci_avg25_incidence <- exp(confint(model_CVD_Incidence_avg25)[, 2])
irr_lci_avg25_incidence
irr_uci_avg25_incidence
model_CVD_Incidence_avg975 <- pglm(CVD_Incidence ~ average_975,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Incidence_avg975)
irr_avg975_incidence <- exp(coef(model_CVD_Incidence_avg975))
irr_avg975_incidence
irr_lci_avg975_incidence <- exp(confint(model_CVD_Incidence_avg975)[, 1])
irr_uci_avg975_incidence <- exp(confint(model_CVD_Incidence_avg975)[, 2])
irr_lci_avg975_incidence
irr_uci_avg975_incidence
model_CVD_Incidence_var25 <- pglm(CVD_Incidence ~ variation_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Incidence_var25)
irr_var25_incidence <- exp(coef(model_CVD_Incidence_var25))
irr_var25_incidence
irr_lci_var25_incidence <- exp(confint(model_CVD_Incidence_var25)[, 1])
irr_uci_var25_incidence <- exp(confint(model_CVD_Incidence_var25)[, 2])
irr_lci_var25_incidence
irr_uci_var25_incidence
### CVD_Prevalence
model_CVD_Prevalence_avg25 <- pglm(CVD_Prevalence ~ average_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Prevalence_avg25)
irr_avg25_prevalence <- exp(coef(model_CVD_Prevalence_avg25))
irr_avg25_prevalence
irr_lci_avg25_prevalence <- exp(confint(model_CVD_Prevalence_avg25)[, 1])
irr_uci_avg25_prevalence <- exp(confint(model_CVD_Prevalence_avg25)[, 2])
irr_lci_avg25_prevalence
irr_uci_avg25_prevalence
model_CVD_Prevalence_avg975 <- pglm(CVD_Prevalence ~ average_975,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Prevalence_avg975)
irr_avg975_prevalence <- exp(coef(model_CVD_Prevalence_avg975))
irr_avg975_prevalence
irr_lci_avg975_prevalence <- exp(confint(model_CVD_Prevalence_avg975)[, 1])
irr_uci_avg975_prevalence <- exp(confint(model_CVD_Prevalence_avg975)[, 2])
irr_lci_avg975_prevalence
irr_uci_avg975_prevalence
model_CVD_Prevalence_var25 <- pglm(CVD_Prevalence ~ variation_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_Prevalence_var25)
irr_var25_prevalence <- exp(coef(model_CVD_Prevalence_var25))
irr_var25_prevalence
irr_lci_var25_prevalence <- exp(confint(model_CVD_Prevalence_var25)[, 1])
irr_uci_var25_prevalence <- exp(confint(model_CVD_Prevalence_var25)[, 2])
irr_lci_var25_prevalence
irr_uci_var25_prevalenceIRR - cvd prevalence incidence mortality
| Humidity | Percentile | Mortality [95% CIs] | Prevalence[95% CIs] | Incidence[95% CIs] |
|---|---|---|---|---|
| Average | 2.5 | 1.01 [1.004,1.017] | 1.016[1.015,1.017] | 1.013[1.01,1.016] |
| Average | 97.5 | 1.002 [1.001,1.003] | 1.002[1.002,1.002] | 1.002[1.001,1.002] |
| Variation | 2.5 | 1.014 [0.996,1.032] | 1.014[1.01,1.017] | 1.011[1.001,1.021] |
PAF
Click here to view/hide the code
gbd_data_uncorrected <- read_csv(here::here("data","gbd_data_uncorrected.csv"))
gbd_data_uncorrected <- gbd_data_uncorrected |>
filter(Year %in% 2011:2019)
# rename state.name as State
#combined_data_temperature_forpanel
#<- combined_data_temperature_forpanel |>
#dplyr::rename(State = state_name,
# Year = year)
# merge temperature data and gbd data
combined_data_uncorrectedgbd <- merge(
combined_data_relhum_forpanel,
gbd_data_uncorrected,
by = c("State", "Year")
)
paf_data_formapviz <- combined_data_uncorrectedgbd |>
dplyr::select(State, Year, CVD_Death, avg_975,
CVD_Prevalence, CVD_Incidence) |>
filter(Year==2019) |>
mutate(paf= round((CVD_Death * 0.002)),
paf_prev= round(CVD_Prevalence * 0.002),
paf_incid=round(CVD_Incidence * 0.002)) |>
dplyr::select(-CVD_Death,-avg_975) |>
#rename State as state
dplyr::rename(state = State)Click here to view/hide the code
library(sf)
library(tidyverse)
india <- readRDS(here::here("data",
"spatial_files",
"India_states.rds"))
india <- india |>
rename(state = NAME_1) |>
# recode orissa as odisha
mutate(state = recode(state, "Orissa" = "Odisha")) |>
# recode uttaranchal as uttarakhand
mutate(state = recode(state, "Uttaranchal" = "Uttarakhand"))
# combine paf_data_formapviz with india based on state variable
india_paf <- india |>
left_join(paf_data_formapviz, by = "state")
india_paf <- india_paf |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# mortality map
india_paf |>
ggplot() +
geom_sf(aes(fill = paf)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf), size = 3, color = "white") +
labs(fill = "CVD deaths",
title = "CVD deaths attributable to high average humidity in Indian states",
subtitle = "Year 2019") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
# Incidence map
india_paf |>
ggplot() +
geom_sf(aes(fill = paf_incid)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_incid), size = 3, color = "white") +
labs(fill = "CVD Incidence",
title = "CVD Incidence attributable to high average humidity in Indian states",
subtitle = "Year 2019") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
# Prevalence map
india_paf |>
ggplot() +
geom_sf(aes(fill = paf_prev)) +
scale_fill_gradient(low = "orange", high = "midnightblue",
labels=scales::comma) +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_prev), size = 3, color = "white") +
labs(fill = "CVD Prevalence",
title = "CVD Prevalence attributable to high average humidity in Indian states",
subtitle = "Year 2019") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
Click here to view/hide the code
## daly
model_CVD_DALY_avg25 <- pglm(CVD_DALY ~ average_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_DALY_avg25)
irr_avg25_DALY<- exp(coef(model_CVD_DALY_avg25))
irr_avg25_DALY
irr_lci_avg25_DALY <- exp(confint(model_CVD_DALY_avg25)[, 1])
irr_uci_avg25_DALY <- exp(confint(model_CVD_DALY_avg25)[, 2])
irr_lci_avg25_DALY
irr_uci_avg25_DALY
model_CVD_DALY_avg975 <- pglm(CVD_DALY ~ average_975,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_DALY_avg975)
irr_avg975_DALY <- exp(coef(model_CVD_DALY_avg975))
irr_avg975_DALY
irr_lci_avg975_DALY <- exp(confint(model_CVD_DALY_avg975)[, 1])
irr_uci_avg975_DALY <- exp(confint(model_CVD_DALY_avg975)[, 2])
irr_lci_avg975_DALY
irr_uci_avg975_DALY
model_CVD_DALY_var25 <- pglm(CVD_DALY ~ variation_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_DALY_var25)
irr_var25_DALY <- exp(coef(model_CVD_DALY_var25))
irr_var25_DALY
irr_lci_var25_DALY <- exp(confint(model_CVD_DALY_var25)[, 1])
irr_uci_var25_DALY <- exp(confint(model_CVD_DALY_var25)[, 2])
irr_lci_var25_DALY
irr_uci_var25_DALY
## yld
model_CVD_YLD_avg25 <- pglm(CVD_YLD ~ average_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_YLD_avg25)
irr_avg25_YLD <- exp(coef(model_CVD_YLD_avg25))
irr_avg25_YLD
irr_lci_avg25_YLD <- exp(confint(model_CVD_YLD_avg25)[, 1])
irr_uci_avg25_YLD <- exp(confint(model_CVD_YLD_avg25)[, 2])
irr_lci_avg25_YLD
irr_uci_avg25_YLD
model_CVD_YLD_avg975 <- pglm(CVD_YLD ~ average_975,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_YLD_avg975)
irr_avg975_YLD <- exp(coef(model_CVD_YLD_avg975))
irr_avg975_YLD
irr_lci_avg975_YLD <- exp(confint(model_CVD_YLD_avg975)[, 1])
irr_uci_avg975_YLD <- exp(confint(model_CVD_YLD_avg975)[, 2])
irr_lci_avg975_YLD
irr_uci_avg975_YLD
model_CVD_YLD_var25 <- pglm(CVD_YLD ~ variation_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_YLD_var25)
irr_var25_YLD <- exp(coef(model_CVD_YLD_var25))
irr_var25_YLD
irr_lci_var25_YLD <- exp(confint(model_CVD_YLD_var25)[, 1])
irr_uci_var25_YLD <- exp(confint(model_CVD_YLD_var25)[, 2])
irr_lci_var25_YLD
irr_uci_var25_YLD
## YLL
model_CVD_YLL_avg25 <- pglm(CVD_YLL ~ average_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_YLL_avg25)
irr_avg25_YLL <- exp(coef(model_CVD_YLL_avg25))
irr_avg25_YLL
irr_lci_avg25_YLL <- exp(confint(model_CVD_YLL_avg25)[, 1])
irr_uci_avg25_YLL <- exp(confint(model_CVD_YLL_avg25)[, 2])
irr_lci_avg25_YLL
irr_uci_avg25_YLL
model_CVD_YLL_avg975 <- pglm(CVD_YLL ~ average_975,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_YLL_avg975)
irr_avg975_YLL <- exp(coef(model_CVD_YLL_avg975))
irr_avg975_YLL
irr_lci_avg975_YLL <- exp(confint(model_CVD_YLL_avg975)[, 1])
irr_uci_avg975_YLL <- exp(confint(model_CVD_YLL_avg975)[, 2])
irr_lci_avg975_YLL
irr_uci_avg975_YLL
model_CVD_YLL_var25 <- pglm(CVD_YLL ~ variation_25,
data = combined_data,
family = poisson,
model = "within",
index = c("State", "Year"))
summary(model_CVD_YLL_var25)
irr_var25_YLL <- exp(coef(model_CVD_YLL_var25))
irr_var25_YLL
irr_lci_var25_YLL <- exp(confint(model_CVD_YLL_var25)[, 1])
irr_uci_var25_YLL <- exp(confint(model_CVD_YLL_var25)[, 2])
irr_lci_var25_YLL
irr_uci_var25_YLLIRR - cvd DALY YLD YLL
| Humidity | Percentile | DALY [95% CIs] | YLD[95% CIs] | YLL[95% CIs] |
|---|---|---|---|---|
| Average | 2.5 | 1.001 [1,1.003] | 1.018[1.012,1.023] | 1.001[0.999,1.002] |
| Average | 97.5 | 1.001 [1.001,1.001] | 1.002[1.001,1.003] | 1.001[1.001,1.001] |
| Variation | 2.5 | 1.004 [1,1.007] | 1.015[0.999,1.03] | 1.003[0.999,1.007] |
PAF for DALY, YLL, YLD
Click here to view/hide the code
gbd_data_uncorrected_economy <- read_csv(here::here("data","CVD_DALY_YLL_YLD_2019_Statewise_Uncorrected.csv"))
# rename state.name as State
#combined_data_temperature_forpanel
#<- combined_data_temperature_forpanel |>
#dplyr::rename(State = state_name,
# Year = year)
# merge temperature data and gbd data
gbd_data_uncorrected_economy_forpaf <- merge(
combined_data_relhum_forpanel,
gbd_data_uncorrected_economy,
by = c("State", "Year")
)
paf_economydata_formapviz <- gbd_data_uncorrected_economy_forpaf |>
dplyr::select(State, Year, avg_975,
DALY,YLL,YLD) |>
filter(Year==2019) |>
mutate(paf_daly= round((DALY * 0.000999001)),
paf_yld= round(YLL * 0.002),
paf_yll=round(YLD * 0.002)) |>
#dplyr::select(-CVD_Death,-avg_975) |>
#rename State as state
dplyr::rename(state = State)Click here to view/hide the code
library(sf)
library(tidyverse)
india <- readRDS(here::here("data",
"spatial_files",
"India_states.rds"))
india <- india |>
rename(state = NAME_1) |>
# recode orissa as odisha
mutate(state = recode(state, "Orissa" = "Odisha")) |>
# recode uttaranchal as uttarakhand
mutate(state = recode(state, "Uttaranchal" = "Uttarakhand"))
# combine paf_data_formapviz with india based on state variable
india_paf_economy <- india |>
left_join(paf_economydata_formapviz, by = "state")
india_paf_economy <- india_paf_economy |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# DALY map
india_paf_economy |>
ggplot() +
geom_sf(aes(fill = paf_daly)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_daly), size = 3, color = "white") +
labs(fill = "CVD DALY",
title = "CVD DALY attributable to high average humidity in Indian states",
subtitle = "Year 2019") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
# yld map
india_paf_economy |>
ggplot() +
geom_sf(aes(fill = paf_yld)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_yld), size = 3, color = "white") +
labs(fill = "CVD YLD",
title = "CVD YLD attributable to high average humidity in Indian states",
subtitle = "Year 2019") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
# YLL map
india_paf_economy |>
ggplot() +
geom_sf(aes(fill = paf_yll)) +
scale_fill_gradient(low = "orange", high = "midnightblue",
labels=scales::comma) +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_yll), size = 3, color = "white") +
labs(fill = "CVD YLL",
title = "CVD YLL attributable to high average humidity in Indian states",
subtitle = "Year 2019") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
Projection maps for 97.5th percentile of average relative humidity
CVD Deaths
Click here to view/hide the code
all_forecasts_deaths <- read_csv(here::here("output", "cvd_deaths","combined_forecasts_all_states.csv"))
all_forecasts_2030_deaths <- all_forecasts_deaths |>
dplyr::filter(Year == 2030) |>
mutate(paf= round((Forecasted_CVD_Deaths * 0.002),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_deaths)Click here to view/hide the code
india <- readRDS(here::here("data",
"spatial_files",
"India_states.rds"))
india <- india |>
rename(state = NAME_1) |>
# recode orissa as odisha
mutate(state = recode(state, "Orissa" = "Odisha")) |>
# recode uttaranchal as uttarakhand
mutate(state = recode(state, "Uttaranchal" = "Uttarakhand"))
# combine paf_data_formapviz with india based on state variable
india_paf_deaths <- india |>
left_join(all_forecasts_2030_deaths, by = "state")
india_paf_deaths <- india_paf_deaths |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# mortality map
india_paf_deaths |>
ggplot() +
geom_sf(aes(fill = paf)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf), size = 3, color = "white") +
labs(fill = "CVD deaths",
title = "CVD deaths attributable to high average relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
CVD Prevalence
Click here to view/hide the code
all_forecasts_prev <- read_csv(here::here("output", "cvd_prev","combined_forecasts_all_states.csv"))
all_forecasts_2030_prev <- all_forecasts_prev |>
dplyr::filter(Year == 2030) |>
mutate(paf= round((Forecasted_CVD_Prevalences * 0.002),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_prev)Click here to view/hide the code
# combine paf_data_formapviz with india based on state variable
india_paf_prev <- india |>
left_join(all_forecasts_2030_prev, by = "state")
india_paf_prev <- india_paf_prev |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# prevalence map
india_paf_prev |>
ggplot() +
geom_sf(aes(fill = paf)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf), size = 3, color = "white") +
labs(fill = "CVD prevalence",
title = "CVD prevalence attributable to high average relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
CVD Incidence
Click here to view/hide the code
all_forecasts_incidence <- read_csv(here::here("output", "cvd_inci","combined_forecasts_all_states.csv"))
all_forecasts_2030_incidence <- all_forecasts_incidence |>
dplyr::filter(Year == 2030) |>
mutate(paf= round((Forecasted_CVD_Incidences * 0.002),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_incidence)Click here to view/hide the code
# combine paf_data_formapviz with india based on state variable
india_paf_inci <- india |>
left_join(all_forecasts_2030_incidence, by = "state")
india_paf_inci <- india_paf_inci |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# incidence map
india_paf_inci |>
ggplot() +
geom_sf(aes(fill = paf)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf), size = 3, color = "white") +
labs(fill = "CVD incidence",
title = "CVD incidence attributable to high average relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
Projection maps for 2.5th percentile of average relative humidity
CVD Deaths
Click here to view/hide the code
all_forecasts_deaths <- read_csv(here::here("output", "cvd_deaths","combined_forecasts_all_states.csv"))
all_forecasts_2030_deaths <- all_forecasts_deaths |>
dplyr::filter(Year == 2030) |>
mutate(paf_deaths= round((Forecasted_CVD_Deaths * 0.0099),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_deaths)Click here to view/hide the code
india <- readRDS(here::here("data",
"spatial_files",
"India_states.rds"))
india <- india |>
rename(state = NAME_1) |>
# recode orissa as odisha
mutate(state = recode(state, "Orissa" = "Odisha")) |>
# recode uttaranchal as uttarakhand
mutate(state = recode(state, "Uttaranchal" = "Uttarakhand"))
# combine paf_data_formapviz with india based on state variable
india_paf_deaths <- india |>
left_join(all_forecasts_2030_deaths, by = "state")
india_paf_deaths <- india_paf_deaths |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# mortality map
india_paf_deaths |>
ggplot() +
geom_sf(aes(fill = paf_deaths)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_deaths), size = 3, color = "white") +
labs(fill = "CVD deaths",
title = "CVD deaths attributable to low average relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
CVD Prevalence
Click here to view/hide the code
all_forecasts_prev <- read_csv(here::here("output", "cvd_prev","combined_forecasts_all_states.csv"))
all_forecasts_2030_prev <- all_forecasts_prev |>
dplyr::filter(Year == 2030) |>
mutate(paf_prev= round((Forecasted_CVD_Prevalences * 0.0157),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_prev)Click here to view/hide the code
# combine paf_data_formapviz with india based on state variable
india_paf_prev <- india |>
left_join(all_forecasts_2030_prev, by = "state")
india_paf_prev <- india_paf_prev |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# prevalence map
india_paf_prev |>
ggplot() +
geom_sf(aes(fill = paf_prev)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_prev), size = 3, color = "white") +
labs(fill = "CVD prevalence",
title = "CVD prevalence attributable to low average relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
CVD Incidence
Click here to view/hide the code
all_forecasts_incidence <- read_csv(here::here("output", "cvd_inci","combined_forecasts_all_states.csv"))
all_forecasts_2030_incidence <- all_forecasts_incidence |>
dplyr::filter(Year == 2030) |>
mutate(paf_inci= round((Forecasted_CVD_Incidences * 0.0128),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_incidence)Click here to view/hide the code
# combine paf_data_formapviz with india based on state variable
india_paf_inci <- india |>
left_join(all_forecasts_2030_incidence, by = "state")
india_paf_inci <- india_paf_inci |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# incidence map
india_paf_inci |>
ggplot() +
geom_sf(aes(fill = paf_inci)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_inci), size = 3, color = "white") +
labs(fill = "CVD incidence",
title = "CVD incidence attributable to low average relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
Projection maps for 2.5th percentile of variation in relative humidity
CVD Deaths
Click here to view/hide the code
all_forecasts_deaths <- read_csv(here::here("output", "cvd_deaths","combined_forecasts_all_states.csv"))
all_forecasts_2030_deaths <- all_forecasts_deaths |>
dplyr::filter(Year == 2030) |>
mutate(paf_deaths= round((Forecasted_CVD_Deaths * 0.0138),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_deaths)Click here to view/hide the code
india <- readRDS(here::here("data",
"spatial_files",
"India_states.rds"))
india <- india |>
rename(state = NAME_1) |>
# recode orissa as odisha
mutate(state = recode(state, "Orissa" = "Odisha")) |>
# recode uttaranchal as uttarakhand
mutate(state = recode(state, "Uttaranchal" = "Uttarakhand"))
# combine paf_data_formapviz with india based on state variable
india_paf_deaths <- india |>
left_join(all_forecasts_2030_deaths, by = "state")
india_paf_deaths <- india_paf_deaths |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# mortality map
india_paf_deaths |>
ggplot() +
geom_sf(aes(fill = paf_deaths)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_deaths), size = 3, color = "white") +
labs(fill = "CVD deaths",
title = "CVD deaths attributable to low variation in relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
CVD Prevalence
Click here to view/hide the code
all_forecasts_prev <- read_csv(here::here("output", "cvd_prev","combined_forecasts_all_states.csv"))
all_forecasts_2030_prev <- all_forecasts_prev |>
dplyr::filter(Year == 2030) |>
mutate(paf_prev= round((Forecasted_CVD_Prevalences * 0.0138),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_prev)Click here to view/hide the code
# combine paf_data_formapviz with india based on state variable
india_paf_prev <- india |>
left_join(all_forecasts_2030_prev, by = "state")
india_paf_prev <- india_paf_prev |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# prevalence map
india_paf_prev |>
ggplot() +
geom_sf(aes(fill = paf_prev)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_prev), size = 3, color = "white") +
labs(fill = "CVD prevalence",
title = "CVD prevalence attributable to low variation in relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())
CVD Incidence
Click here to view/hide the code
all_forecasts_incidence <- read_csv(here::here("output", "cvd_inci","combined_forecasts_all_states.csv"))
all_forecasts_2030_incidence <- all_forecasts_incidence |>
dplyr::filter(Year == 2030) |>
mutate(paf_inci= round((Forecasted_CVD_Incidences * 0.01088),0)) |>
dplyr::rename(state = State)
rmarkdown::paged_table(all_forecasts_2030_incidence)Click here to view/hide the code
# combine paf_data_formapviz with india based on state variable
india_paf_inci <- india |>
left_join(all_forecasts_2030_incidence, by = "state")
india_paf_inci <- india_paf_inci |>
mutate(centroid = st_centroid(geometry)) %>%
mutate(x = st_coordinates(centroid)[,1],
y = st_coordinates(centroid)[,2])
# incidence map
india_paf_inci |>
ggplot() +
geom_sf(aes(fill = paf_inci)) +
scale_fill_gradient(low = "orange", high = "midnightblue") +
# Add paf on the states as text with x and y coordinates
geom_text(aes(x = x, y = y, label = paf_inci), size = 3, color = "white") +
labs(fill = "CVD incidence",
title = "CVD incidence attributable to low variation in relative humidity in Indian states",
subtitle = "Year 2030") +
theme_classic() +
# Make x-axis and y-axis blank
theme(axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.line = element_blank())